---
title: "kindergarten spatial distance analysis at Salzburg, Austria"
output:
  flexdashboard::flex_dashboard:
    storyboard: yes
    social: menu
    source: embed
    css: style.css
---

```{r setup, include=FALSE}
library(ggplot2)
library(plotly)
library(leaflet)
library(htmltools)
library(mapview)
library(sf)

set.seed(10)
# Read the data
salzburg <- read_sf("data/salzburg_city.shp.shp")
kinder <- read_sf("data/kindergarden.shp.shp") %>% 
  cbind(st_coordinates(.$geometry))  
motor <- read_sf("data/motorway.shp.shp")
kinder$PLAETZE_sqrt  <- sqrt(kinder$PLAETZE)*1.5
min_distance <- sapply(seq_len(nrow(kinder)), function(x) min(st_distance(kinder[x,], motor)))
kinder$min_distance <- min_distance
kinder$classes <- kmeans(kinder$min_distance, 4)$cluster

# Ordering classes
fct_class <- as.factor(kinder$classes)
levels(fct_class) <- c("Low", "Very low", "Medium", "High")
fct_class <- factor(fct_class, c("Very low", "Low", "Medium", "High"))
kinder$classes <- fct_class

kinder$pop_message_map01  <- paste(
  sep = "
",
  sprintf("name: %s", kinder$NAME),
  sprintf("seats: %s", kinder$PLAETZE)
)

kinder$pop_message_map02  <- paste(
  sep = "
",
  sprintf("name: %s", kinder$NAME),
  sprintf("distance: %s meters", round(kinder$min_distance))
)

kinder$pop_message_map03  <- paste(
  sep = "
",
  sprintf("name: %s", kinder$NAME),
  sprintf("class: %s", kinder$classes),
  sprintf("seats: %s", kinder$PLAETZE)
)

pal <- colorNumeric(
  palette = "YlOrRd",
  domain = kinder$min_distance
)

pal_class <- colorFactor(
  palette = "RdYlBu",
  domain = kinder$classes,
  reverse = TRUE,
  n = 4
)
```

### **Kindergarten's** spatial distribution at Salzburg, Austria. {data-commentary-width=400}

```{r}
leaflet(kinder) %>% 
  setView(13.03982, 47.79800, 12) %>% 
  addProviderTiles("Esri.WorldImagery", group = "ESRI basemap") %>%
  addProviderTiles(providers$CartoDB.Voyager, group = "Carto basemap") %>%
  addPolylines(data = salzburg, color = "#696969", opacity = 1, weight = 2) %>% 
  addPolylines(data = motor, color = "red", group="Motorways") %>% 
  addCircleMarkers(lng = ~X, 
                   lat = ~Y, 
                   weight = 1.5,
                   radius = ~PLAETZE_sqrt,
                   color = "black",
                   opacity = 1,
                   fillOpacity = 0.1,
                   group = "Kindergartens",
                   popup = ~pop_message_map01) %>% 
  addLayersControl(
    baseGroups = c("Carto basemap", "ESRI basemap"),
    overlayGroups = c("Motorways", "Kindergartens"),
    options = layersControlOptions(collapsed = FALSE)
  )
```

------------------------------------------------------------------------


Kindergartens always need to be in the safest possible areas within a city. However, this is not always possible. 


As you can see on the left map, some kindergartens are exposed to the proximity of motorways. 

Making use of the distance analysis we will determine what **kindergartens are more exposed** to the risk of proximity to highways.

**The radius of circles** is related to the number of seats of kindergartens.



### How close to the **motorways** are they?. {data-commentary-width=400}

```{r}
leaflet(kinder) %>% 
  setView(13.03982, 47.79800, 12) %>% 
  addProviderTiles("Esri.WorldImagery", group = "ESRI basemap") %>%
  addProviderTiles(providers$CartoDB.Voyager, group = "Carto basemap") %>%
  addPolylines(data = salzburg, color = "#696969", opacity = 1, weight = 2) %>% 
  addPolylines(data = motor, color = "red", group="Motorways") %>% 
  addCircleMarkers(lng = ~X, 
                   lat = ~Y, 
                   weight = 1.5,
                   color = "black",
                   fillColor = ~colorQuantile("YlOrRd", min_distance)(min_distance),
                   opacity = 1,
                   fillOpacity = 1,
                   group = "Kindergartens",
                   popup = ~pop_message_map02) %>% 
  addLayersControl(
    baseGroups = c("Carto basemap", "ESRI basemap"),
    overlayGroups = c("Motorways", "Kindergartens"),
    options = layersControlOptions(collapsed = FALSE)
  ) %>% 
  addLegend("bottomleft", pal = pal, values = ~min_distance,
            title = "Distance in meters",
            opacity = 1
  )
```

------------------------------------------------------------------------


Points in **light yellow** represent the kindergartens **closest** to the motorways while the **red** ones are the **furthest away**. See the legend at the bottom left.


### Classification according to the proximity of the **motorways**. {data-commentary-width=400}

```{r}
leaflet(kinder) %>% 
  setView(13.03982, 47.79800, 12) %>% 
  addProviderTiles("Esri.WorldImagery", group = "ESRI basemap") %>%
  addProviderTiles(providers$CartoDB.Voyager, group = "Carto basemap") %>%
  addPolylines(data = salzburg, color = "#696969", opacity = 1, weight = 2) %>% 
  addPolylines(data = motor, color = "red", group="Motorways") %>% 
  addCircleMarkers(lng = ~X, 
                   lat = ~Y, 
                   radius = ~PLAETZE_sqrt,
                   weight = 1.5,
                   color = "black",
                   fillColor = ~colorFactor("RdYlBu", reverse = TRUE, classes)(classes),
                   opacity = 1,
                   fillOpacity = 1,
                   group = "Kindergartens",
                   popup = ~pop_message_map03) %>% 
  addLayersControl(
    baseGroups = c("Carto basemap", "ESRI basemap"),
    overlayGroups = c("Motorways", "Kindergartens"),
    options = layersControlOptions(collapsed = FALSE)
  ) %>% 
  addLegend("bottomleft", pal = pal_class, values = ~classes,
            title = "",
            opacity = 1
  )
```

------------------------------------------------------------------------


Using the _proximity to highways_ values from the map above, the kindergartens have been classified as high, medium, low, and very low risk. 



### Statistic for each class. {data-commentary-width=400}

```{r}
plaetze_classes <- classInt::classIntervals(kinder$PLAETZE, 4, "quantile", dataPrecision = 1)
plaetze_cut <- cut(kinder$PLAETZE, plaetze_classes$brks, include.lowest = TRUE)
plaetze_cut_order <- factor(plaetze_cut, levels(plaetze_cut), ordered = TRUE)
kinder$PLAETZE_class <- plaetze_cut_order

p <- ggplot(data = kinder, aes(x = classes, fill = PLAETZE_class)) +
  geom_bar(position = "dodge") + 
  guides(fill=guide_legend("PLAETZE")) +
  xlab("") +
  ylab("Number of kindergarten")
ggplotly(p)
```

------------------------------------------------------------------------



In total there are **73 kindergartens** in the Salzburg district, of which **14** are **Very Low**, **19 Low**, **21 Medium**, and **19 High** risks

As the graph on the left shows, kindergartens classified as **High risk** also have the **highest number of students**. This pattern could be related to the population growth that occurs, most of the time, near to important highways.



### Kindergardens with at least 25 seats. {data-commentary-width=400}

```{r}
leaflet(kinder[kinder$PLAETZE >25,]) %>% 
  setView(13.03982, 47.79800, 12) %>% 
  addProviderTiles("Esri.WorldImagery", group = "ESRI basemap") %>%
  addProviderTiles(providers$CartoDB.Voyager, group = "Carto basemap") %>%
  addPolylines(data = salzburg, color = "#696969", opacity = 1, weight = 2) %>% 
  addPolylines(data = motor, color = "red", group="Motorways") %>% 
  addCircleMarkers(lng = ~X, 
                   lat = ~Y, 
                   radius = ~PLAETZE_sqrt,
                   weight = 1.5,
                   color = "black",
                   fillColor = ~colorFactor("RdYlBu", reverse = TRUE, classes)(classes),
                   opacity = 1,
                   fillOpacity = 1,
                   group = "Kindergartens",
                   popup = ~pop_message_map03) %>% 
  addLayersControl(
    baseGroups = c("Carto basemap", "ESRI basemap"),
    overlayGroups = c("Motorways", "Kindergartens"),
    options = layersControlOptions(collapsed = FALSE)
  ) %>% 
  addLegend("bottomleft", pal = pal_class, values = ~classes,
            title = "",
            opacity = 1
  )
```

------------------------------------------------------------------------


Considering only the kindergartens with more than 25 seats.  We can observe that the number of kindergartens was reduced by 57% (from 73 to 31). Furthermore, clearly shown that the kindergartens furthest from the cities are the most affected by the filter applied.



### Kindergardens with at least 25 seats. **Area of influence** {data-commentary-width=400}


```{r}
buffer <- st_buffer(endCapStyle="ROUND",kinder[kinder$PLAETZE>=25,], 750/111111, nQuadSegs = 10) %>% 
  st_union() %>% 
  st_intersection(salzburg)

leaflet(kinder[kinder$PLAETZE >= 25,]) %>% 
  setView(13.03982, 47.79800, 12) %>% 
  addProviderTiles("Esri.WorldImagery", group = "ESRI basemap") %>%
  addProviderTiles(providers$CartoDB.Voyager, group = "Carto basemap") %>%
  addPolylines(data = salzburg, color = "#696969", opacity = 1, weight = 2) %>% 
  addPolylines(data = salzburg, color = "#696969", opacity = 1, weight = 2) %>% 
  addPolylines(data = motor, color = "red", group="Motorways") %>%
  addPolygons(data = buffer, 
              color = "black",
              fillColor = "#134ff2", 
              weight = 1,
              group = "Area of influence") %>% 
  addCircleMarkers(lng = ~X, 
                   lat = ~Y, 
                   radius = 4,
                   weight = 1.5,
                   color = "black",
                   opacity = 1,
                   fillOpacity = 1,
                   group = "Kindergartens") %>% 
  addLayersControl(
    baseGroups = c("Carto basemap", "ESRI basemap"),
    overlayGroups = c("Motorways", "Kindergartens", "Area of influence"),
    options = layersControlOptions(collapsed = FALSE)
  )
```

------------------------------------------------------------------------


Taking 750 meters as a constant value for the range of influence, it can be seen that the kindergartens' services area covers 50% of the Salzburg district.

Although there are better ways to estimate the radius of influence than the euclidean distance (e.g. cost distance or network distance), this information is still very useful for decision makers to know **where to build kindergartens**. 

To find suitable locations, we also recommend taking into account population density and proximity to parks.